'''
Created on Nov 6, 2009

@author: mkiyer
'''
import sys
import collections
import os
import glob
import logging
import numpy as np
import operator
from parse_cufflinks_results import parse_transcripts_tmap, parse_transcripts_expr
from stats import calc_copa_score
from cancer_specificity import count_cancer_outliers, get_library_params

if __name__ == '__main__':
    logging.basicConfig(level=logging.DEBUG)
    results_path = sys.argv[1]
    tmap_files = glob.glob(os.path.join(results_path, "*.tmap"))    
    exprs = collections.defaultdict(lambda: {})
    libraries = []
    library_params = get_library_params()    

#    stop_after = 2
    for i, tmap_file in enumerate(tmap_files):
#        if i == stop_after:
#            break        
        logging.debug('processing %s' % tmap_file)
        library_id = os.path.splitext(os.path.basename(tmap_file))[0]
        libraries.append(library_id)        
        for tx in parse_transcripts_tmap(open(tmap_file)):
            if (tx.class_code == '=') or (tx.class_code == 'c'):
                if tx.major_iso_id == tx.cuff_id:
                    library_covs = exprs[tx.ref_id]
                    library_covs[library_id] = tx.rpkm
    
#    outfhd = open('ref_transcripts.txt', 'w')
#    outfhd.write('ref_id')
#    outfhd.write('\trecurrences\tcopa_75\tcopa_90\tcopa_99')
#    for library in libraries:
#        outfhd.write('\t%s' % library)
#    outfhd.write('\n')
    outfhd = sys.stdout
    
    recurrences = collections.defaultdict(lambda: 0)
    for ref_id, library_covs in exprs.iteritems():        
        # get list of library ids and coverages for this reference transcript
        library_ids = library_covs.keys()
        covs = [library_covs[x] for x in library_ids]
        # compute statistics for this locus
        recurrences[len(covs)] += 1    
        score75 = calc_copa_score(np.array(covs), r=0.75)
        score90 = calc_copa_score(np.array(covs), r=0.90)
        score99 = calc_copa_score(np.array(covs), r=0.99)
        noutliers, ncancer = count_cancer_outliers(library_ids, covs, library_params)
        if noutliers == 0:
            frac_cancer = 0
        else:
            frac_cancer = float(ncancer) / noutliers

        if noutliers <= 1:
            median = covs[-1]
        else:
            median = np.median(sorted(covs)[len(covs)-noutliers:])
          

        # print the results
        line = '\t'.join([str(len(covs)),
                          '%.3f' % (frac_cancer),
                          str(noutliers),
                          str(ncancer),
                          '%.3f' % (score75),
                          '%.3f' % (score90),
                          '%.3f' % (score99),
                          '%s' % (ref_id),
                          '%.3f' % median,
                          ','.join(library_ids),
                          ','.join(map(str,covs))])
        outfhd.write('%s\n' % line)
        
#        outfhd.write('%s\t%d\t%.3f\t%.3f\t%.3f' % (ref_id, len(library_covs.values()), 
#                                                   score75, score90, score99))
#        for library_id in libraries:
#            if library_id in library_covs:
#                outfhd.write('\t%.3f' % library_covs[library_id])
#            else:
#                outfhd.write('\tNone')
#        outfhd.write('\n')
    logging.debug('recurrences: %s' % (recurrences))    

 
    


